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ABSTRACT 


The  nature  of  generated  machine  error  in  finite  digital  calculations 
is  discussed.      The  arithmetic  of  range  numbers  is  developed,    and 
examples  are  given  demonstrating  the  use  of  range  arithmetic  as  a 
tool  for  automatic  error  analysis.     A  computer  program  is  developed, 
utilizing  the  TYPE  OTHER  feature  of  FORTRAN-63  in  conjunction  with 
the  CDC- 1 604  digital  computer,    which  enables  the  user  to  perform 
automatic  error  analysis  during  computation,    and  a  number  of  programs 
are  presented  using  this  feature. 
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1.       Introduction. 

The  utilization  of  high-speed  digital  computers  to  perform  lengthy 
and  complex  mathematical  computations  is  widespread.     Unfortunately, 
the  user  often  assumes  that  the  computer  output  will  necessarily  be 
accurate  if  his  input  data  is  accurate  and  his  program  logic  is  sound. 
This  need  not  be  the  case,    and  such  an  assumption  may  well  lead  the 
user  to  incorrect  or  highly  inaccurate  results.      The  difficulty  arises 
because  of  the  very  nature  of  digital  computation. 

A  digital  computer  has  a  finite  memory,    and  cannot  carry  within  it 
infinite  precision  numbers,    or  perform  infinitely  accurate  calculations. 
Each  number  stored  in  a  digital  computer  must,    of  necessity,    contain 
a  finite  number  of  digits.     As  a  result,    during  each  step  of  computation, 
the  computer  automatically  truncates  or  rounds  off  numbers  to  conform 
with  its  memory  storage  requirements.     Each  arithmetic  operation 
introduces  an  error  in  precision  due  to  such  truncation,    which  may 
accumulate  with  other  such  errors  to  produce  sizeable  inaccuracies  in 
results.      It  is  also  possible  that  the  errors  may  cancel  and  give  a  fairly 
accurate  result;  but,    in  either  case,    the  user  will  have  no  knowledge  of 
whether  or  not  any  error  had  accumulated  in  computation. 

It  would  be  desirable  to  eliminate  precision  errors  entirely,    but 
hardware  limitations  make  it  impossible  to  do  so.      It  is,    however, 
possible  to  estimate  and  control  the  errors  arising  in  finite  digital 
computations,    by  appropriate  choice  of  an  algorithm.      One  method  is 


called  range  arithmetic,    or  interval  arithmetic    [1,  2  J.      This  method 
enables  the  generated  error  to  be  analysed  automatically  within  the 
computer,    during  computation. 

The  technique  of  range  arithmetic  replaces  any  real  number  of  infi- 
nite precision  by  two  numbers  of  finite  precision,    one  of  which  is  a 
lowerbound  of  the  real  number,    and  the  other  is  a  corresponding  upper - 
bound  of  the  real  number.      We  can  then  be  certain  that  the  real  number 
lies  in  the  range  of  the  two  finite  precision  numbers. 

In  this  thesis,    we  will  discuss  the  arithmetic  of  range  numbers,    the 
application  of  range  arithmetic  to  digital  computation,    and  present  a 
computer  program  designed  for  use  on  the  CDC-1604  computer  in  con- 
junction with  FORTRAN-63,    which  will  provide  a  means  of  automatic 
error  analysis  in  digital  computation. 
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2.       The  Arithmetic  of  Range  Numbers. 

In  this  section,    we  will  discuss  the  properties  of  range  numbers, 
and  the  arithmetic  operations  pertaining  to  them. 

Let  x  be  any  real  number.     Associated  with  this  real  number  is  a 
range  number  X,    which  has  the  following  properties: 

(1)  X  is  a  closed  bounded  interval; 

(2)  X  contains  x. 

Then  define  the  range  number  X,    associated  with  the  real  number  x, 
to  be  the  closed  bounded  interval 

X   =    [xL,  xU], 
where  xL  and  xU  are  defined  by  the  transformations  L  and  U  as  follows: 

L 

X >  XL 


x — -— >xU 


We  may  then  write 

L,  U 


•>  X   =    [xL,  xU]       . 


The  transformations  L  and  U  are  determined  by  the  precision  allowed 
in  computation.      If,    for  example,    we  would  like  to  find  the  range  number 
associated  with   e  =  2.  71828 ,    and  we  are  allowed  only  two  sig- 
nificant figures,    then 

X   =    [2.7,  2.8]    . 
On  the  other  hand,    if  allowed  five  significant  figures,    then 

X   =   [2.  7182  ,    2.7183]  . 
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For  integers  (known  with  exact  precision),    the  range  number  will  have 
the  same  upper  and  lower  bound,    equal  to  the  integer  itself.      Then,    if 
x   =    3,     X   =    [3,    3]  . 

We  can  see,    then,    that  the  range  number  X  will  represent  the 
smallest  closed  interval  containing  x,    for  a  given  set  of  transformations 
L  and  U. 

The  arithmetic  operations  associated  with  range  numbers  can  be 
defined  in  the  following  manner.      Let  (op)  represent  an  arithmetic 
operation,    and  consider  two  range  numbers 

X    =    [xL,    xU]    , 
and 

Y    =   [yL,    yU]    . 
Then, 

Z    =   X  (op)    Y      =   [zL,    zU]    , 
where 

zL   =   min  [xL  (op)  yL  ,    xL(op)yU,    xU(op)yL,    xU(op)yU] 
zU   =   max[xL  (op)  yL  ,    xL(op)yU,    xU(op)yL,    xU(op)yU]    . 
It  can  be  seen  that  the  real  number    z    =    x  (op)  y   will  be  in  the  interval  Z, 
and  that  Z  represents  the  smallest  interval  that  contains  z. 
Now  let 

(  +  )     -     denote  range  addition 
(-)     -     denote  range  subtraction 
(*)     -     denote  range  multiplication 
(/)     -     denote  range  division 
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From  the  definition  above,    it  can  be  seen  that 

( 1 )  Z    =    X  (+)    Y    =   [  xL  +  yL  ,    xU  +  yU  ] 

(2)  -X    =    [-xU,    -xL] 

(3)  Z    =    X(-)    Y   =   X(  +  )(-Y)   =   [xL-yU,   xU  -  yL  ] 

(4)  X    =     Y    if  and  only  if 
xL    =    yL  ,    and  xU   =   yU 

(5)  Z    =    X  (*)  Y,    where 

zL  =  min  [  xL  *  y L  ,  xL  *  yU  ,  xU  *  yL  ,  xU  *  yU  ] 
zU  =  max  [  xL  *  yL  ,  xL  *  yU  ,  xU  *  yL  ,  xU  *  yU  ] 

The  values  of  zL  and  zU  can  also  be  obtained  using  a  table  of  sign 

dis  criminations . 

Using  the  table  on  page  14,    we  need  only  calculate  one  product  for 

each  end  point,    except  in  sign  discrimination  5. 

(6)  Z    =    X  (/)  Y,    where 

zL  =  min  [  xL  /  yL  ,  xL  /  yU  ,  xU  /  yL  ,  xU  /  yU  ] 
zU  =  max  [  xL  /  yL  ,  xL  /  yU  ,  xU  /  yL  ,  xU  /  yU  ] 

Note,    that  in  the  case  in  which  Y  brackets  zero,    i.  e.  ,    Y  =  [  -  1  ,    3], 
range  division  is  undefined.      This  is  equivalent  to  forbidden  division  by 
zero  in  real  arithmetic. 
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We  now  will  give  some  illustrative  examples  of  computations  using 
range  arithmetic. 

[a,    b]   (  +  )  [c,    c]    =    [a  +  c,   b  +  c] 

[a,    b]   (*)  [c,    c]    =    [  ca  ,    cb  ]  for  c  £  0 

[a,    b]   (*)  [  c,    c]    =    [  cb  ,    ca  ]  for  c  £  0 

[a,    b]   (-)  [a,    b]    =    [a-b,    b-a] 

[a,    b]  (/)  [a,    b]    =    [a/b,    b/a]  for  a  >  0 

[1,    2]  (*)  [-4,    -3]    =    [-8,    -3] 
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3.       Range  Arithmetic  and  Automatic  Error  Analysis. 

In  finite  digital  computation,    errors  are  generated  in  the  process 
of  round-off,    truncation,    and  normalization.      One  has  no  method  of 
determining  the  seriousness  of  such  errors  unless  some  sort  of  error 
analysis  is  performed.     A  method  is  necessary  that  can  evaluate  the 
error  generated  at  each  step  in  the  computation,    and  also  keep  track  of 
the  accumulated  error.      The  techniques  of  range  arithmetic  are  suited 
to  this  purpose.     If,    in  addition  to  doing  real  arithmetic,    we  parallel  the 
computation  with  range  arithmetic,    the  range  numbers  calculated  will 
keep  track  of  the  error,    and  the  final  result  will  be  the  smallest  closed 
interval  containing  the  "true",    infinite  precision  result. 

Error  analysis  via  range  arithmetic  not  only  furnishes  information 
as  to  the  size  of  machine-generated  errors,    but  it  also  can  be  utilized 
to  evaluate  the  accuracy  of  alternative  computational  schemes.      Indeed, 
the  error  generated  in  one  type  of  computation  may  be  quite  different 
from  that  obtained  using  another  method  of  computation.     A  few  examples 
will  illustrate  this  point. 

Example  1 . 

Consider  the  problem  of  finding  the  roots  of  the  equation 

2 
x      +    100, OOOx  +    100    =    0    . 

We  will  concern  ourselves  with  finding  the  larger  root  of  this 

equation,    and  we  will  assume  that  the  computer  can  carry  six  significant 
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figures.     We  will  first  calculate  the  result  using  ordinary  arithmetic 
and  the  quadratic  formula,    and  then  use  range  arithmetic  in  conjunction 
with  the  quadratic  formula.     The  notation  will  be  that  which  is  equivalent 
to  normalized  floating-point  operations  on  a  computer.     For  example, 
the  number  100,  000  will  be  denoted  by  .  100000E6,    where  E6  denotes  the 
power  of  ten  multiplying  .  100000. 
The  equation 


2 

ax     +   bx  +   c   =   0 


has,    as  its  larger  root, 


♦  J*z 


x    =    (  -b  +    A/b     -  4ac)  /  2a    . 
We  then  calculate  our  solution: 


-.100000E6+    a/.100000E11    -    .400000E3 
X    "  . 2E1 

-.100000E6  +    //.  100000E11    -    .000000E11 

.  2E1 

-.100000E6+    .  100000E6 


.  2E1 
000000E00 


or 


x    =    0 

Now,    in  order  to  do  the  same  problem  in  range  arithmetic,   we  assume 
some  uncertainty  in  the  coefficients  of  our  original  equation.      We  then 
replace  the  original  coefficients  by  range  numbers  reflecting  the  un- 
certainty.    For  example,   we  will  replace  100,000  by  [99999.9,  100,001] 

17 


The  equation  to  be  solved  now  is 

[.999999E0,    .  100001E1  ]  x2    +    [.999999E5,    .100001E6]x 


Now, 


b2  = 


ac  = 


4ac  = 


b      -  4ac  = 


/ 


b     -  4ac  = 


-b  \N  b      -  4ac 


and 


x  = 


+    [.999999E2,    .100001E3]     =    0   . 

.999999E5,    .100001E6]2 
.  99999E10,    .  100003E11] 

.999999E0,    .  100001E1]  x  [  .  999999E2,    .100001E3] 
. 999998E2,     .100003E3] 

.  399999E3,    .400012E3] 

.  999998E10,  .  100003E11]  -  [.  399999E3,    .  400012E3] 
.999997E10, . 100003E11] 

. 999998E5,    .100002E6] 

.999998E5,    .100002E6]    -[.999999E5,    .100001E6] 
-.120000E4,     .210000E4] 


J1-.120000E4,    .21000QE4] 
Q399999E3,    .  400012E3J 


=    £.299994E1,    .  525002ElJ 

By  using  range  arithmetic,    we  have  kept  accurate  error  bounds  and  have 
found  that  our  root  is  in  the  interval   [  -  2.  9,    5.  3],    or 
[  -  .  299994E1,     .  525002E1]  .      Our  result  tells  us  that  the  true  value 
could  be  anywhere  within  the  calculated  interval.      However,    the  interval 
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is  wide,    and  could  contain  the  value  zero,   which  was  our  answer  using 
real  arithmetic.      We  would  like  to  tighten  the  interval,    so  that  we  may 
be  more  confident  of  the  computer  results.      The  next  example  will 
illustrate  how,    with  a  little  knowledge  of  sources  of  error  generated  in 
finite  calculations,    we  can  obtain  more  satisfactory  results. 

Example  2. 

Consider  the  formula  for  the  larger  root  of  a  quadratic  equation, 

x  =   (-b+A/b     -  4ac)  /  2a 
and  divide  numerator  and  denominator  by  b  so  that 

x   =    (  -  1  +     JT~-~$)  b  /  2a 

where 

s   =   4ac /b 

Expanding  the  radical  in  a  power  series,    we  obtain 

2,_  3 


J\   -  s  =  1   -  s/2  -  s    /8  -  s    /16  - 
The  equation  for  x  becomes 


x  =  -  (s/2  +  s2/8  +  s3/l6  + )b/2a 


Now,    for  the  equation 


2 
x     +  100,  OOOx  +  100  =  0 


we  get 

(.  400000E1)  (.  100000E1)  (.  100000E3) 


s      = 


(.  100000E6)2 


=    .400000E-7 
s/2     =    .  200000E-7 
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Then 


or, 


s    /8     =     .  200000E-15 


2 
b    /2a     =    (.  100000E6)  /  (.  200000E1)  (.  100000E1) 


=  . 500000E5 

x     =  -(.  500000E5)  (.  200000E-7  +  .  200000E-15) 

=  -(.  500000E5)  (.  200000E-7) 

=  -    . 100000E-2 

x     =  -    . 001 


We  now  solve  the  same  problem  using  range  arithmetic,    with  the 
modified  quadratic  formula. 

Our  equation  is 

[.999999E0,    .  100001E1]  x2  +   [.999999E5,    .100001E6]x 

+     [.999999E2,    .  100001E3]     =    0   . 
In  example  2,    we  found  that 

b       =    [.999998E10,     .  100003E11] 

4ac     =    [.399999E3,       .400012E3] 

therefore, 

[.399999E3,        .400012E3] 
S     ~    [.999998E10,     . 100003E11] 

=    [.399987E-7,     .400013E-7] 

s/2     =    [.399987E-7,     .  400013E-7]  /  [  .  200000E1 ,     .  200000E1] 
=    [.199994E-7,     .  200007E-7] 
20 


and 


2 
s        =    [.  159987E-14,     .  160009E-14] 

2 
s    /8     =    [.  199983E-15,     .200012E-15] 

2 
s/2  +  s    /8     =    [.199994E-7,     .  200008E-7] 

b/a     =    [.999999E5,     .  100001E6]/ [.  999999E0,     .  100001E1] 
=    [.999989E5,     .  100002E6] 

b/2a     =    [.499994E5,     .  500005E5] 

x     =    -  [.499994E5,    .  500005E5]  [  .  199994E-7,    .  200008E-7] 
=    [-.  100005E-2,     -.999958E-3] 

x     =    [-.00100005,      -.000999958]    . 

Thus,    we  can  see  that  the  method  in  example  2  generates  far  less 
error  than  the  straight  application  of  the  quadratic  formula.     Indeed, 
the  previous  examples  illustrate  the  use  of  range  arithmetic  in  evaluating 
a  given  computational  technique. 


or, 
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4.       QRANGE7    -    Program  Description  and  Usage. 

The  TYPE  OTHER  capability  in  the  FORTRAN-63  compiler  for  the 
CDC- 1604  computer  allows  the  relatively  simple  implementation  of 
range  arithmetic  for  the  user.     Using  the  QRANGE7  package,    the 
programmer  need  only  declare  all  floating-point  variables  to  be 
TYPE  RANGE7(3),    and  follow  the  usual  rules  of  FORTRAN  programming. 
The  QRANGE7  package  supplies  the  necessary  subroutines  for  range 
arithmetic  computations. 

Range  arithmetic  calculations  are  made  carrying  a  triple  of  numbers 
(xF,    xL,    xU),    rather  than  a  double  as  stated  previously.      The  first 
number  in  the  triple  is  the  ordinary  floating-point  result,    as  calculated 
using  standard  floating-point  arithmetic.      The  remaining  numbers  in  the 
triple  are  the  lower  and  upper  range  numbers,    respectively.     Note  that 
xL  ^  xF  ^  xU  .      In  any  range  calculation,    the  "real"  part  of  the  triple 
is  computed  last  and  remains  in  the  accumulator,    so  that  transfers  and 
comparisons    (such  as  the  IF  statement)   will  follow  the  same  logical 
branches  as  in  ordinary  floating-point  arithmetic. 

Standard  floating-point  arithmetic  cannot  be  used  on  either  lower 
or  upper  range  numbers  since,    in  floating-point  calculations,    round-off 
and  normalization  are  done  automatically,    and  would  reduce  the  accu- 
racy of  the  range  interval.      In  QRANGE7,    the  arithmetic  performed  on 
range  numbers  is  un-normalized,    and  separate  subroutines  are  provided 
which  truncate  the  lower  range  numbers  and  round  up  by  one  bit  in  the 
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least  significant  position,    or  truncate,    the  upper  range  numbers  as 
required.      The  un-normalized  arithmetic  is  performed  by  unpacking  the 
operands,    putting  them  in  un-normalized  form,    doing  the  appropriate 
arithmetic,    re -normalizing  and  repacking  the  result. 

Input/output  may  be  accomplished  in  one  of  two  ways.     The 
programmer  may  use  EQUIVALENCE  statements  to  transfer  input  /  output 
of  range  numbers  as  real  numbers,    where  each  range  number  is  equiva- 
lent to  three  real  numbers,    or  he  may  use  the  Input  /  Output  for 
Multi-Word  TYPE  OTHER  package  [3]  and  transfer  the  range  numbers 
directly  in  and  out  of  the  computer  memory.      The  following  examples 
will  illustrate  the  use  of  these  methods,    and  the  use  of  QRANGE7  . 


A. 


TYPE  RANGE7(3)    A7,  B7,  X7 

DIMENSION    A7(3)  ,  AR(9)  ,  B7(3)  ,  BR(9)  ,  X7(3,  3)  ,  XR(27) 

EQUIVALENCE    (A7,  AR)  ,  (B7,  BR)  ,  (X7,  XR) 

READ    10,  AR,  BR 
10    FORMAT(8F10.0) 

DO    20    1=1,  3 

DO    20    J=l,  3 
20    X7(J,  I)  =  A7(I)*B7(J) 

DO    30    1=1,  25,  3 

30    PRINT    40,XR(I)  ,XR(I  +  1)  ,XR(I+2)  ,XR(I)  ,  XR(I  +  1),     XR(I+2) 

40    FORMAT(//15X,  3(E20.  10,  10X)/15X,  3(020,  10X)) 

END 
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In  this  example,    we  have  assumed  that  the  programmer  reads  in  the 
three  parts  of  the  range  numbers  using  real  floating-point  numbers  that 
are  simultaneously  stored  as  the  three  parts  of  the  range  numbers  by 
use  of  the  EQUIVALENCE  statement.      The  next  example  will  illustrate 
the  same  program,    but  using  the  INPUT  /  OUTPUT  for  MULTI-WORD 
TYPE  OTHER  package   [3], 


B. 


TYPE    RANGE7(3)    A7,  B7,  X7 

DIMENSION    A7(3)  ,  B7(3)  ,X7(3,  3) 

CALL    READ7 

READ    10,A7,B7 
10    FORMAT(8F10.0) 

DO    20    1=1,3 

DO    20    Jcl,  3 

X7(J,  I)    =   A7(I)*B7(J) 

CALL    PRINT7 
20    PRINT    30,    X7(J,  I) 
30    FORMAT(//15X,  3(E20.  10,  10X)  /  15X,  3(020,  10X)) 

END 
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5.       Matrix  Inversion  Using  Range  Arithmetic. 

In  order  to  test  QRANGE7  ,    and  evaluate  the  errors  generated  in 
a  particular  program,   a  number  of  matrix  inversions  were  run.     In  one 
case,    the  matrix  to  be  inverted  was  generated  randomly,    using  the 
RANF(-l)  library  function  of  FORTRAN-63  .     About  15  programs  were 
run  in  this  manner.     In  the  second  case,   the  input  matrix  was  specified, 
and  the  individual  elements  of  the  matrix  were  assumed  to  be  inaccurate. 
Thus,    they  were  inputed  as  range  numbers.     Both  cases  used  the 
MATINV2  subroutine  provided  by  the  U.  S.    Naval  Postgraduate  School 
computer  facility,    in  conjunction  with  QRANGE7  .      Two  range  sub- 
routines had  to  be  used:    ABS7,    which  takes  the  absolute  value  of  a 
range  number;    and  Q0Q06700,   which  complements  the  range  accumulator 
[5].     (A  listing  of  these  subroutines  appears  in  Appendix  II.)     The 
coding  for  both  cases  was  as  follows: 

Case  1.        Randomly  Generated  Matrices 
PROGRAM  RANMAT 
TYPE    RANGE7(3)    XMAT  ,  XINV 
DIMENSION    XMAT(5,  5)  ,  XINV(5,  5)  ,  XJUNK(2000) 
READ    2000  ,  N 
DO    90    1=1  ,  N 
90      XJUNK(I)    =   RANF(-l) 
DO    10    1=1,5 
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DO    10    J=l,5 

B(I,  J)    =    RANF(  -  1)*1000 
10      XMAT(I,  J)    =    B(I,  J) 

CALL    PRINT7 

DO    20    J=l,  5 

DO    20    1  =  1,  5 
20      PRINT    2001  ,  XMAT(I,  J)  ,  XMAT(I,  J) 

CALL    MATINV2(XMAT  ,  XINV  ,  5,  5) 

CALL    PRINT7 

DO    30    J=l,5 

DO    30    1  =  1,5 
30      PRINT    2001  ,XINV(I,  J),  XINV(I,  J) 

2000  FORMAT(F10.  0) 

2001  FORMAT(  /  / 1  5X  ,  3(E20.  10  ,  10X)  /  15X  ,  3(020,  10X)) 
END 

Case  2.        Input  Matrix  Containing  Range  Elements 
PROGRAM    MMATRIX 
TYPE    RANGE7(3)    A,  B 
DIMENSION    A(3,  3)  ,B(3,  3) 
CALL  READ7 
DO    10    J=l,3 
DO    10     1  =  1,3 
10      READ    2000,  A( I,  J) 
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CALL    PRINT7 

DO    20    J=l,3 

DO  20  1  =  1,3 
20      PRINT    2001  ,  A(I,  J),  A(I,  J) 

CALL    MATINV2(A,  B,  3,  3) 

CALL  PRINT 7 

DO    30    J=l,3 

DO  30  1  =  1,3 
30      PRINT    2001  ,B(I,  J),  B(I,  J) 

2000  FORMAT(3F10.  0) 

2001  FORMAT(//15X  ,  3(E20.  10,  10X)  /  15X,  3(020,  10X) ) 
END 

The  results  obtained  using  random  matrices  were  highly  accurate. 
For  5x5  matrices,    seven  decimal  places  of  accuracy  were  obtained, 
and  the  range  of  error  was  small.     It  might  be  noted  that  the  conversion 
of  the  real  random  matrix  to  a  matrix  of  range  numbers  introduced  no 
errors,    and  hence  the  resulting  range  intervals  on  the  inverse  elements 
were  solely  a  consequence  of  digital  computation.      On  the  other  hand, 
the  error  generated  in  the  second  case  was  quite  large,    and  decimal 
accuracy  varied  from  zero  to  one.      In  fact,   large  errors  were  generated 
for  input  coefficients  with  ranges  of  approximately  1.  OE-2.      It  was  found, 
in  both  cases,    that  the  error  generated  increased  with  increasing 
matrix  dimension.      This  is  to  be  expected,    since  the  number  of 
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calculations  done  in  inverting  the  matrix  is  increased.     For  example,    a 
12x12  random  matrix  was  inverted,    and  it  was  found  that  only  four  to 
five  significant  figures  could  be  expected,    vice  the  seven  decimal  accu- 
racy in  the  5x5  case. 

One  my  circumvent  the  inaccuracies  obtained  when  the  input 
elements  are  range  numbers,    by  the  following  device  [4].     Suppose  we 
wish  to  invert  X,    whose  elements  are  range  numbers.      Take  X     to  be 
the  matrix  of  mid-points  of  the  elements  of  X,    so  that  X  =  X     +  E  . 
The  elements  of  E  are  range  numbers  of  the  form  [  -e  ,    e]  .     Obtain 
X  using  range  arithmetic,    and  let 

E'     =    [I-X0-1E    +(X0_1E)2     -    ]    X0_1    (-EXq"1). 

If  the  elements  of  E  are  small  compared  to  X  ,    then  the  series  will 

converge,    and  one  can  write 

X"  l     =    X      +    E;  . 

The  results  of  a  typical  run  for  an  inversion  of  a  randomly  generated 
5x5  matrix  are  shown  on  the  following  page. 
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Element 

Real  Part 

L 

ower 

Upper 

(1.  1) 

4.  6396189281E01 

Same 

as  Real 

Sam< 

i  as  Real 

(1.2) 

8. 6916184500E02 

(1.3) 

3.  6039163981E02 

(1.4) 

5.  8772778565E02 

(1,  5) 

9.4693749440E02 

(2,  1) 

5.  5609401546E02 

(2,2) 

8.9089112203E02 

(2,3) 

4. 0143081011E02 

(2,4) 

4. 2098029158E02 

(2,5) 

6.  1093177836E02 

(3,  1) 

9.9636585643E02 

(3,2) 

1.  6340010236E02 

(3,  3) 

4.  665803721 5E02 

(3,4) 

5.  0479887773E02 

(3,  5) 

2.  0507283101E02 

(4,  1) 

2.  7500284988E02 

(4,2) 

4.  8510492363E02 

(4,3) 

2.4488146363E02 

(4,4) 

4.  1884968178E02 

(4,  5) 

1.9965179004E02 

(5,  1) 

8.  7792113349E02 

(5,2) 

2. 3254672355E02 

(5,3) 

3. 5002222285E00 

1 

(5,4) 

3. 2092384145E02 

(5,5) 

6.4308479314E02 

I 

\ 

\ 

/ 

X  -  MATRIX 


TABLE  2 
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We  will  present  the  real  part,    and  the  least  significant  digits  of  the 
lower  and  upper  range  numbers  of  X  -  INVERSE.     The  exponents  will  not 
be  repeated. 


Element 

Real  Part 

Lower 

Upper 

(1,  1) 

-9.  1422040         728E-04 

854 

60  7 

(1,2) 

7.419425         7112E-04 

6975 

7249 

(1,3) 

9.3418389         435E-05 

070 

864 

(1.4) 

-5.779019         5190E-05 

5709 

4688 

(1.5) 

6.2948717         498E-04 

477 

517 

(2,1) 

-1.122627           041E-03 

052 

032 

(2,2) 

1.6441862         954E-03 

944 

965 

(2,3) 

-1.3201931          162E-03 

165 

159 

(2,4) 

1.4714936         870E-03 

865 

875 

(2,5) 

5.5237176         832E-05 

720 

995 

(3,1) 

1.0620289         855E-03 

838 

873 

(3,2) 

1.1383145         460E-03 

437 

478 

(3,3) 

2.3018511          807E-03 

800 

812 

(3,4) 

-3.6993050        965E-03 

973 

956 

(3,5) 

-2.2307817         966E-03 

966 

962 

(4,  1) 

6.4802494         170E-04 

040 

338 

(4,2) 

-3.0018154         973E-03 

990 

955 

(4,3) 

-5.102677         3200E-05 

3669 

2679 

(4,4) 

4.0549999         820E-03 

811 

825 

(4,5) 

6.  5487719         463E-04 

433 

477 

(5,1) 

1.3248530         755E-03 

742 

766 

(5,2) 

-1.1561055         861E-04 

968 

706 

(5,3) 

3.6280024         837E-04 

794 

865 

(5,4) 

-2.4566800         524E-03 

531 

520 

(5,5) 

3.6100559         502E-04 

481 

523 

X  -  INVERSE 

FIGURE  3 
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6.       Conclusions. 

Range  arithmetic  is  a  powerful  tool,    not  only  for  error  estimation 
in  computations  of  a  scientific  or  engineering  nature,    but  for  the  evalu- 
ation and  comparison  of  alternative  numerical  algorithms  as  well.     The 
QRANGE7   package  provides  the  user  with  a  simple  means  for  such 
analysis.      The  drawback  in  the  use  of  QRANGE7  ,    at  this  point,    is  that 
computing  time  is  increased.      There  are  many  points  where  the  program 
could  be  increased  in  efficiency,    notably  in  the  Q1Q04770    and   Q1Q05770 
subroutines  for  range  multiplication  and  range  division.      The  author 
feels,    however,    that  even  with  the  increase  in  computing  time,    the  pay- 
off in  the  use  of  range  arithmetic  is  so  great  that  it  more  than  outweighs 
the  disadvantage. 

A  few  function  subprograms  have  been  provided  for  use  in  conjunc- 
tion with  QRANGE7  .     They  are  used  exactly  as  one  would  use  the  library 
subroutines  in  FORTRAN,    except  a  "7"  replaces  the  "F"  at  the  end  of 
the  function  name.      One  must  remember  to  declare  these  functions 
TYPE    RANGE7(3)   before  execution.     Much  work  remains  to  be  done  in 
this  area,    and  it  is  hoped  that,    eventually,    a  complete  library  of  func- 
tions will  be  available  for  use  with  QRANGE7  . 

It  might  be  noted,    at  this  point,    that  there  are  errors  generated  in 
the  conversion  of  decimal-to-binary  and  binary-to-decimal  numbers 
during  input/output,    and  these  have  not  been  taken  into  account.      The 
author  feels  that  these  errors  are  negligible  compared  to  the  machine 
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generated  error  occurring  during  computation. 

The  results  obtained  using  QRANGE7  indicated  that  the  technique  of 
range  arithmetic  can,    indeed,    be  used  to  provide  error  information  at 
each  stage  of  computation.      Unfortunately,    the  author  did  not  have  the 
time  to  try  it  out  on  some  of  the  more  common  numerical  algorithms, 

such  as  the  Runge  -  Kutta  method  for  solving  differential  equations,    or 

th 
the  Newton  -  Raphson  method  for  determining  roots  of  n      order 

polynomials.     It  is  hoped  that  this  will  be  done  in  the  future. 

In  conclusion,    it  is  felt  that  range  arithmetic  is  of  such  value  that 

it  warrants  inclusion  in  future  modifications  of  FORTRAN  or  other 

algebraic  compilers  as  a  standard  TYPE. 
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APPENDIX  I 
QRANGE7   PROGRAM  LISTING  [  5] 

Entry  Points: 

Q1Q00770,  Q1Q01770,  Q1Q02770,  Q1Q03770, 
Q1Q04770,  Q1Q05770,  Q1Q10770,  Q1Q00710, 
Q1Q01710,  Q1Q02710,  Q1Q03710,  Q1Q04710, 
Q1Q05710,  Q1Q10710,  Q1Q10170,  Q1Q00700, 
Q1Q01700,  Q1Q02700,  Q1Q03700,  Q1Q04700, 
Q1Q05700,  Q1Q10700,  Q1Q10070,  UNPACK7, 
RNGAD7,  RNGMU7,  RNGDI7,  QUIRKL7,  QUIRKU7, 
REPACK7,    ERROR777 
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APPENDIX  II 
SUBROUTINES  FOR  USE  WITH   QRANGE7    [5] 


1.  ABS7 

2.  SQRT7 

3.  INT7 

4.  Q0Q06700 

5.  Q2Q07770 
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